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SECTION  I 


INTRODUCTION 

The  major  uncertainty  in  the  calculation  of  the  conductivity  at  a 
given  temperature  in  certain  plasmas  representative  of  those  produced  as 
combustion  products  in  MHD  generators  is  due  to  the  uncertainty  in  the 
low  energy  momentum  transfer  cross  sections  for  electron  scattering  in 
water  vapor. ^  The  most  reliable  values  of  such  low  energy  cross  sections 
are  those  deduced  from  electron  swarm  experiments.  The  determination  of 
cross  sections  from  swarm  experiments  involves  solving  the  Boltzmann  equation 
for  the  electron  energy  distribution  function  using  some  trial  set  of 
assumed  cross  sections.  Certain  integrals  of  the  distribution  function  yield 
the  quantities  measured  in  swarm  experiments,  such  as  drift  velocities  and 
diffusion  coefficients.  By  comparing  the  calculated  values  of  these 
measurable  quantities  with  those  found  experimentally,  it  is  possible  to 
reach  useful  conclusions  about  the  trial  set  of  cross  sections.  An  iterative 
procedure  can  then  allow  deduction  of  a  set  of  cross  sections  which  is  consistent 
with  swarm  observations.  The  low  energy  cross  sections  are  thus  determined 
and  are  some  of  the  most  reliable. 

The  conventional  method  of  solving  the  Boltzmann  equation  in  the  context 
of  electron  swarm  experiments  relies  on  the  approximation  that  the  inelastic 
scattering  cross  sections  are  much  smaller  than  the  elastic.  Because  of  the 
large  dipole  moment  in  H^O,  the  rotational  cross  sections  will  be  rather  large 
and  comparable  to  the  elastic.  Therefore,  before  attempting  to  determine 
cross  sections  in  H^O  from  swarm  experiments  it  has  been  necessary  to  improve 
on  the  solution  and  to  evaluate  the  errors  under  some  circumstances. 

This  report  details  our  efforts  in  that  direction  and  shows  results  for 
several  cases  using  the  extended  analyses. 
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The  central  problem  Is  to  deduce  the  distribution  function  from 
assumptions  about  the  cross  sections.  The  solution  must  be  done  numerically 
and  in  the  course  of  the  analysis  it  must  be  done  many  times.  It  is  necessary 
to  strike  a  balance  between  accuracy  and  required  computer  resources. 

The  numerical  problem  is  sufficiently  difficult  that  it  is  necessary  to 
Involve  assumptions  which  are  of  questionable  validity  such  as  the  assumption 
of  a  nearly  spherical  electron  energy  distribution  function.  A  major  concern 
is  the  evaluation  of  the  errors  generated.  Some  of  the  errors  arise  because 
of  physical  assumptions  which  have  the  effect  of  changing  the  problem. 

In  evaluating  the  error  it  is  necessary  to  judge  whether  the  simplified 
problem  is  similar  enough  to  the  Initial  problem.  Another  broad  class 
of  errors  arises  from  inadequacies  in  the  numerical  analysis.  The  usual 
situation  is  that  these  errors  depend  on  both  the  procedure  and  the  data. 

For  both  classes  of  error  we  have  made  substantial  progress  in  working 
out  procedures  with  less  error  than  before. 
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SECTION  II 


THE  BOLTZMANN  EQUATION 


The  motion  of  electrons  traveling  through  a  neutral  gas  under  the  influence 
of  an  electric  field  may  be  described  mathematically  by  the  Boltzman  equation, 

|£(r,v,t)  +  v*$rf (r,v,t)  +  a  •  V^f(r,v,t)  =  c[f(r,v,t)]  (1) 

where  f(r,v,t)  is  the  electron  energy  distribution  function.  The  acceleration, 
a,  is  that  produced  by  the  field  (a  =  eE/m) .  The  right  hand  side  of  the  equation 
is  the  rate  of  change  in  the  distribution  function  due  to  collisions.  Since 
the  only  external  force  in  the  problem  is  a  uniform  electric  field,  the  distri¬ 
bution  function  is  cylindrically  symmetric  about  the  field  direction. 

Normally  it  is  assumed  the  experimental  situation  is  adequately  described 
with  the  assumption  that  f(r,v,t)  is  independent  of  position  and  time.  These 
assumptions  xeaa  to  tne  nrst  two  terms  ot  Eq.  (1)  being  zero.  Actually,  the 
configuration  of  the  experiments  imply  spatial  gradients  which  can  also 
cause  temporal  variations.  The  impact  of  the  density  gradients  on  the 

distribution  function  is  relatively  small  and  can  be  found  by  a  perturbative 

2 

scheme  proposed  by  Skullerud.  First,  we  discuss  the  spatially  independent 
equation  and  its  solutions. 

With  the  assumption  that  the  first  two  terms  of  Eq.  (1)  are  negligible 


the  Boltzmann  equation  becomes 


f(0) 


(v)  -  c[fw'(v)]  =  o  . 


where  the  collision  term  has  the  form. 


c[f(v)]  -  -NvQ(v,6  )f(v)  +  4  1  [  v'f(v’)Q.  (v\6e)6(v'-g(v))dV  .  (3) 

8  v  k-o J  K  3 
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The  first  term  on  the  right  of  Eq.  (3)  represents  scattering-out  of  the  velocity 
3 

volume  d  v  element  centered  about  v  and  the  second  term  is  the  scattering-in 

3 

component  from  all  the  volume  elements  d  v1  centered  about  v'.  The  symbols 

have  the  following  meanings; 

N  *  neutral  particle  density 

Q(v,6  )  *  total  scattering  cross  section 
s 

6g  «  scattering  angle 

Q,(v,6  )  *  scattering  cross  section  for  the  process, 

xC  S 

k=0  for  elastic  scattering  . 

The  delta  functions  serve  the  purpose  of  picking  out  those  electrons  of 

3  t  h 

velocities  vf  that  are  scattered  into  d^v  after  a  collision  of  the  k  type. 

The  integral  over  vf  may  be  immediately  performed  given  g^(v),  the  initial  velocity 
v1  in  terms  of  the  final  velocity  v.  The  g^(v)*s  will  depend  on  the  particular 
processes  considered. 


Since  the  energy  loss  in  elastic  collisions  is  small,  the  second  terra  with 

2 

k=0  may  be  expanded  in  a  Taylor  series  about  v  and  terms  beyond  0(m/M)  may  be 
neglected.  The  resulting  collision  term  is 


*Jl 


Nv  If(o)f,v)  -  f(oo,v)]Q  (v,6  idco* 

o  s 


+n%s£!(i- 


cos0  )Q  (v,0  )v  ffo’jVjdu)1 
s  o  s 

.\2 


(4) 


+  E  Nvj[f  (oj*  (v'  »6S)  “  f  (<»,v)Qk^  ^(v,6s)]dio'  . 


The  sum  is  to  be  made  over  all  inelastic  and  superelastic  processes  with  cross 

fc 

sections  and  Q  ,  respectively. 

In  order  to  convert  the  vector  dependence  of  the  single  function  f  into  the 

more  easily  handled  scalar  dependence  of  a  coupled  set  of  functions  f^ 

we  expand  f  in  terms  of  Legendre  functions  and  the  unknown  coefficients  f^ 

n 

f(v)  -  E  f. (v)P. (cosG)  .  (5) 

Jl-o  *  * 
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Fot  zero  fields  the  distribution  function  must  be  spherically  symmetric. 

At  higher  fields  more  terms  are  needed  to  represent  the  anisotropy  in  f (v) . 

Most  work  in  the  field  has  been  done  assuming  that  n=l  yields  an  adequate 
solution.  This  is  referred  to  later  as  the  two  term  solution.  We  do  not  constrain 
n  to  be  small  but  in  general  find  that  n  need  not  be  larger  than  7. 

After  substituting  this  expansion  in  Eq.  (2) ,  we  multiply  both  sides  by 
P^(cos6)  and  integrate  over  all  6  to  yield  n+1  c  coupled  differential 
equations  for  j=0,  ...  n.  The  coupling  is  between  any  three  adjacent 
expansion  coefficients  according  to  the  equation. 


f  df,  i(e)  ,1  1  i+i  r  dfi+i 

-  E  2j-l  If  de  '2  fj-l(e)J  +  2j+3  ^  de 


(e) 


P.  (cos6)Q  (e,6  )do)' 
.1  S 

(*) 


-NEQ(e)fj(e)  +  Nef  j  (e) 

+N  l  (e±ek)f j(£l?k)  Pj(cos6g)Qk'‘  ''(e±ek.6s)d(jj' 
-  recoil  =  0  j=0,  ...  n 


(6) 


2 

In  writing  this  equation  we  have  made  a  change  of  variables  to  £  =  1/2  mv 
and  have  made  use  of  the  properties  of  Legendre  functions.  These  equations 
are  quite  general  so  far  and  both  superelastic  and  inelastic  as  well  as 
anisotropic  scattering  may  be  explicitly  included.  In  going  from  Eq.  (2) 
to  Eq.  (6)  a  partial  differential  equation  in  two  variables  has  been  converted 
to  Np  ordinary  differential  equations  in  one  variable.  The  principal  merit 
to  this  step  comes  from  the  fact  that  an  adequate  solution  can  be  found 


with  Np  relatively  small. 


SECTION  III 


METHOD  OF  SOLUTION 


The  set  of  first  order  differential  equations  (6)  can  be  compactly  expressed 
in  matrix  form  as 

Al(e)Ff  (e)  +  AO(e)F(e)  +  E  AR(e+e  k)F(e+ek)  »  0  (7) 

where  the  elements  of  the  column  vector  F(e)  are  the  functions  f^e)  and  the 

Afs  are  square  matrices  of  order  N  .  The  A's  are  matrices  whose  rows  contain 

P 

the  coefficients  appearing  in  Eq.  (6).  A1  and  AO  are  tridiagonal  and  the  Ak's 
are  purely  diagonal. 

The  "non-local11  terms,  terms  written  as  functions  of  appearing 

in  Eqs.  (6)  and  (7)  are  due  to  inelastic  and  superelastic  collisions.  It  is  the 
presence  of  these  terms  which  preclude  the  use  of  common  stepping  methods  of 
solution  for  differential  equations.  An  approach  (backward  prolongation^  taken 

3 

by  Frost  and  Phelps,  was  to  begin  by  assuming  a  solution  for  high  energies, 
usually  taken  to  be  the  solution  for  purely  elastic  scattering,  and  then  stepping 
down  to  lower  energies.  If  superelastic  collisions  are  important  it  is  necessary 
to  know  the  solution  at  low  energies  in  order  to  calculate  it  at  high  energies. 

A  global  approach  is  required  in  order  to  include  the  effect  of  the  "non-local" 

terms  in  a  general  way.  The  particular  method  we  have  adopted  is  a  type  of 

4 

collocation  called  Galerkin’s  method. 

To  implement  Galerkin’s  method,  we  make  the  additional  expansion 


N 

F(e)  = 

2  C.S  (e) 
j=l  J  J 

(8a) 

with  F(£)  = 

F(£)  +  6(e)  . 

(8b) 
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C.  is  a  vector  of  length  N  and  each  element  is  by  definition  independent 
J  P 

of  energy*  The  Sj(e)’s  appearing  in  the  expansion  may  be  members  of  any 

convenient  set  of  functions  which  is  complete  or  nearly  complete  in  the  energy 

range  of  interest.  6(e)  is  a  small  error  accounting  for  the  fact  that  a 

finite  number  of  S^’s  do  not  form  a  complete  set  of  functions.  In  our  work 

we  choose  the  S^'s  to  be  'cubic  B-splines.^  Cubic  B-splines  are  third  degree 

polynomials  with  continuous  first  and  second  derivatives.  A  basis  of  such 

functions  is  local  in  the  sense  that  at  every  energy  point  only  a  fixed  number 

(  four  in  the  case  of  cubic  splines)  of  B-splines  is  non-zero.  B-splines 

are  evaluated  quite  easily  from  their  definition  as  a  divided  difference  of 

the  truncated  power  function^  on  some  arbitrary  energy  grid,  . 

s 

A  consequence  of  this  choice  of  basis  functions  is  that  the  distribution 

function  is  zero  above  .  It  is  known  that  the  distribution  is  small  at  high 

Ns 

energies,  and  we  can  always  choose  an  energy  such  that  above  it  the  distribution 
function  is  so  small  we  are  not  interested.  An  important  point  to  examine  is 
that  the  calculated  distribution  function  is  independent  of  the  high  energy 
cut  off  point. 

By  combining  (8a)  and  (7), 


B 

l  Cj [Al(e)Sj (e)  +  A0(e)S  (e)  +  £  A^eie^S  (e±ek)] 
J**  k 


A(e) 


(9) 


where 


A(e)  ■  -[A1(e)6'(e)  +  A0(e)6(e)  +  £  A^eIe^A  (e±£k)  ] 

k 


If  A(e)  were  zero  for  all  e,  then  the  C^'s  could  be  determined  and  would 
yield  a  solution  T(e)  ■  F(e).  It  is  possible  to  approximate  the  ideal  case 
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-- 


of  vanishing  A (c)  by  forcing  A(e)  to  be  outside  the  space  spanned  by  { S ^ } , 
that  is  by  making  A(e)  orthogonal  to  each  Sj.  Therefore  we  require 

<S1(c)|A(c)>  -  0  for  all  i-l,...,N  .  (10) 


Using  relation  10  in  Eq.  (9),  we  obtain 


s 

l  c  [<S  (e)|Al(e)|S!(e)>  +  <S. (e) |A0(e) |s. (e)> 
j-1  11  1  1  1 

+  <S1(c)|i  A^(e±e^) |Sj (e±e^)>]  =  0 


(ID 


or  in  matrix  form, 

H  C  *  0  (12) 

with  equal  to  the  bracketed  term  in  Eq.  (11).  The  matrix  M  is  a  square 

matrix  of  order  N  x  N  but  is  rather  sparse  because  only  three  splines 
p  s 

are  non-zero  at  any  energy  e. 

Eq.  (10)  is  a  set  of  homogeneous  algebraic  equations.  To  make  the  solution 
definite  a  further  statement  is  required  involving  normalization.  It  is 
satisfactory  to  initially  set 

fo<e=°)  =  1  •  (13) 

After  solving  the  matrix  equation  the  distribution  function  is  renormalized 
to  the  condition 

j  e1/2  fo(e)dc  =  1  . 

The  advantages  of  this  method  of  solution  are  many.  Most  importantly, 
the  non-local  terms  appearing  in  Eq.  (11)  and  which  are  due  to  anisotropic 
scattering  for  i>0  are  easily  incorporated  into  the  method.  The  calculation 


of  the  elements  of  the  matrix  M  are  easy  to  perform  using  three  to  nine 
point  Gauss-Legendre  quadrature  techniques.  Many  of  the  calculations  including 
the  method  for  solving  Eq.  (12)  lend  themselves  well  to  the  vector  processing 
options  available  on  the  CRAY-1  computer  that  is  being  used  for  these  calculations. 
The  time  required  is  thereby  reduced  to  seconds  per  run.  This  method  also 
allows  for  an  unevenly  spaced  energy  grid.  However,  since  the  error  increases 
with  the  ratio  of  the  maximum  to  the  minimum  energy  interval  spanned  by 

any  one  B-spline,  some  care  must  be  taken  to  optimize  the  grid  selection. 

1/2 

We  commonly  use  e  spacing.  This  seems  to  work  quite  well,  although  we  suspect 
that  it  is  not  the  optimum  choice. 

The  main  drawback  to  this  method  of  solution  is  in  the  amount  of  computer 
storage  required  for  the  elements  of  the  matrix  M.  No  effort  has  been  made 
so  far  to  take  advantage  of  the  spare ity  of  M  and  the  storage  could  be  reduced 
ny  factors  or  two  to  ten  or  more  wnen  tms  is  aone.  Most  of  the  computing 
time  is  taken  in  determining  the  Al,  AO  and  A^  matrix  elements.  A  significant 
reduction  in  time  could  be  achieved  by  using  some  analytic  form  for  cross 
sections  to  avoid  the  large  numbers  of  interpolations  required  in  the 
case  of  numerical  cross  section  input. 

The  truncation  of  the  Legendre  series  must  be  made  after  an  even  number 
of  terms  if  recoil  energy  losses  are  neglected.  This  is  because  the  number 
degrees  of  freedom  is  equal  to  the  differential  order  of  the  system  only  for  an 
even  number  of  equations  if  there  is  no  f^’  term  (recoil)  in  the  i-0  equation. 

When  recoil  is  included,  the  f0’  term  is  present  but  with  a  very  small  coefficient. 
To  avoid  any  stability  problems  that  may  occur,  we  always  keep  an  even  number 
of  terms  in  the  expansion. 


The  number  of  splines  necessary  to  represent  the  distribution  function 
adequately  depend  on  the  assumed  data.  We  typically  use  about  one  hundred, 
although  with  some  careful  placement  of  the  grid  points,  this  number  may 
be  reduced.  The  choice  of  an  optimum  grid  in  collocation  methods  of  this  sort 
is  a  current  topic  of  research  among  mathematicians  and  there  are  some  mathe¬ 
matical  guidelines  for  choslng  the  grid/*  We  intend  to  implement  some  of  these 
ideas  in  the  computer  code  and  expect  to  cut  down  the  number  of  energy  grid 
points  needed  and  therefore  the  size  of  the  matrix  M. 
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SECTION  IV 


TRANSPORT  COEFFICIENTS 

The  transport  coefficients  are  the  measurable  quantities  in  swarm  experiments. 
It  is  these  parameters  that  we  wish  to  calculate  to  compare  with  the  experiments 
and  hence  determine  the  set  of  cross  sections  that  give  the  best  fit  between 
experiment  and  theory. 

Once  the  distribution  function  f(v)  has  been  calculated,  the  drift  velocity,  W, 
is  simply  the  component  in  the  field  direction  of 

<v>  “  | v  f (v)d^v 

or  <V2>  “  ®  J v cos  0 f (v)d^v 
Because  of  the  orthogonality  of  the  Legendre  functions 

W  *  -j  |  V  f1(v)d3v 

Computation  of  the  diffusion  coefficients  requires  including  the  spatial 

dependence  in  the  distribution  function.  Since  spatial  gradients  of  the 

electron  density  are  small  under  conditions  of  swarm  experiments,  the 

2 

spatial  dependence  may  be  included,  following  Skullerud,  in  a  perturbative  way 
■>  -► 

by  expanding  f(r,v,t)  in  powers  of  the  gradient  of  the  electron  density, 

-  l  f(k)<?)  •  (-?>k  n(r,t) 

k=0  (15) 

-  f^°^(v)n(r,t)  -  ?^(v)  •  Vn(r,t)  +  ... 
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When  this  is  substituted  Into  the  Boltzmann  equation 


n(r,t)  -  C[f*°\v)]j 

+  f(0)(^)  +  {f(0)(^^  +  :.y(1)(v) 

-  Cl?(1)(v)]}  •  $n(r,t) 


+  .  .  . 


The  time  dependence  of  the  electron  density  can  be  expressed  in  terms  of  the 


transport  coefficients,  o)v  , 


=  l  :(k)  •  (-V)kn(r,t) 
3t  k=0 


where 


a)V  ;  =  0  —  no  attachment  or  ionization 
=  W  —  drift  velocity 

(2) 

u>v  '  *=  g  —  diffusion  tensor 

The  first  order  solution  already  discussed  is  obtained  by  assuming  f(r,v,t)  =  f(v), 


Then  from  Eq.  (16), 


a  •  Vvf(0)(v)  -  C[f (0) (v) ]  -  0 


Once  Eq.  (2)  has  been  solved  for  f^(v),  the  second  order  solution  Is  obtained 
by  keeping  two  terms  in  expansion  (15).  The  equation  to  solve  for  t(1)(v)  Is 

«.??(1)(v)  -  C[?(1)(v)]  -  <v-$)f(0)(v)  .  (18) 


■ Si  "it 


wmimtmww 


This  Is  a  vector  equation  and,  because  of  the  cylindrical  symmetry,  there 

ere  only  two  independent  components,  one  for  f|^(v)  and  the  other  for 

f^^(v)  ■  f(1)tf).  Eq.  (18)  can  be  solved  by  the  method  outlined  above 
x  y 

as  the  right  hand  side  is  known. 

The  diffusion  tensor  D  may  be  written  as 

D  -  |  jv?(1)(v)d3v 

The  components  of  interest  are 

D^N  *  [v  f^(v)d^v  —  diffusion  parallel  to  the  field  (19a) 
r  J  z  z 

D^N  *  Jvxf^^(v)d^v  —  diffusion  perpendicular  to  the  field  (19b) 

It  is  D^/U  and  Dp/y  that  are  the  measured  quantities  in  the  experiments 
and  the  relation  to  DTN  or  DpN  is 

£  .  £S%!fi  (20) 

y  w 

In  what  follows,  our  convergence  criterion  is  based  on  values  of  W  and  D^/y 
as  a  function  of  the  number  of  Legendre  functions  in  the  distribution  function 
expansion  rather  than  the  distribution  function  itself  since  they  are  the 
measurables.  We  have  not  yet  spent  much  effort  on  investigating  the  convergence 
of  Dp/y  .  DT/y  is  more  commonly  measured  and  the  convergence  of  Dp/y  is 
expected  to  be  similar  to  DT/y. 

From  the  equations  it  can  be  seen  that  W  and  D^/y  are  functions  of  E/N 
and  the  neutral  gas  composition  only.  (This  is  in  fact  verified  experimentally 
up  to  rather  large  neutral  densities  where  the  Boltzmann  equation  begins 
to  break  down.)  By  varying  E/N,  the  experimenter  is  varying  the  average 
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energy  of  the  electron  swarm.  When  unfolding  data  from  swarm  experiments  to 
determine  cross  sections,  experiments  and  calculations  are  compared  over 
as  wide  a  range  of  E/N  as  possible,  from  near  thermal  average  energies 
to  the  region  of  appreciable  energy  loss  to  Ionization, 
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SECTION  V 


MODEL  RESULTS 

Many  of  the  properties  of  the  solution,  in  particular  the  convergent 
properties,  can  be  seen  most  easily  in  simple  model  cases.  In  this  sectlm,  we 
discuss  the  solutions  for  two  model  situations,  a  model  atom  and  a  methane model, 
each  with  one  Inelastic  cross  section,  before  going  to  the  case  of  a  real  (as, 
nitrogen  (N2>. 

1.  Model  Atom 

The  slmpliest  model  system  we  have  looked  at  is  that  representative  of  a 

model  atom.  The  cross  sections  for  this  model  are  illustrated  in  Fig.  1. 

“16  2 

The  elastic  cross  section  is  constant  at  6*  10  cm  and  the  inelastic 

cross  section  is  of  the  form  k(e-e^),  a  ramp  with  a  threshold  energy  at  e^. 

The  scattering  is  purely  isotropic  and  superelastic  collisions  are  not  induced. 

The  neutral  mass  was  taken  to  be  4  a.m.u.'s.  The  distribution  function,  W 

and  D^/u  were  calculated  for  this  model  over  a  range  of  E/N  from  1  to  48  Td 
-17  2 

(10  V  cm  ),  thresholds  at  .2  and  2  eV,  and  slopes  for  the  inelastic 

—16  2 

cross  section  of  2  and  10  x 10  cm  /eV.  These  values  of  input  parameters 
are  characteristic  of  distribution  functions  with  average  energies  from  a  small 
fraction  to  several  times  the  threshold  energy. 

Since  the  calculated  values  of  transport  coefficients  are  the  quantities 
compared  with  experiment  to  determine  the  cross  sections,  it  is  these 
quantities  with  which  we  are  most  concerned.  In  particular,  we  would  like  to 
evaluate  the  error  introduced  by  truncating  the  Legendre  expansion  for  the 
distribution  function  after  two  terms  as  is  usually  done.  To  this  end 
we  have  calculated  values  of  W  and  D^N  using  as  many  terms  in  the  Legendre  expansion 
as  are  required  for  convergence. 


CROSS  ! 


Figure  1.  Model  atom  cross  sections 
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The  results  of  the  calculations  of  W  and  D^/u  as  a  function  of  E/N  are 

-lfc  2 

shown  for  Np»2  and  Npe6  in  Fig.  2  for  a  value  of  k  *  10x10  cm  /eV.  The  two-term 
values  are  always  greater  than  for  Np*6  and  the  difference  between  D^/p  for  Np=2 
to  Np*4  was  always  much  greater  than  the  changes  from  Np=4  to  the  largest  value  of 
Np  used.  This  indicates  that  the  values  of  the  transport  coefficients  have 
converged  as  a  function  of  Np. 

We  have  compared  values  of  W  and  DT/p  with  Monte  Carlo  calculations. 

The  comparison  with  the  six-term  results  is  within  the  quoted  accuracy^ 

of  the  Monte  Carlo  results  over  the  entire  range  of  parameters  investigated. 

Although  there  are  no  physical  significances  to  the  Legendre  expansion 

coefficients  beyond  the  first,  it  is  instructive  to  compare  their  magnitudes. 

Figure  3  is  a  plot  of  the  f^'s  for  i«0  to  5  at  an  E/N  of  24Td.,  =  *2  eV,  and 

~16  2 

k  *  10  x 10  cm  eV.  This  is  a  case  where  the  two  term  approximation  clearly 

creates  aown.  The  interesting  features  here  are  the  zeros  in  the  higher  order 
coefficients  occurring  at  successively  higher  energies  and  the  maxima  separated 
by  approximately  the  same  energy  for  1^2.  The  most  prominent  feature  is  the 
similarity  in  the  shapes  of  the  f^'s  on  the  high  energy  side  of  the  zeros  of  the 
functions.  The  changes  in  f^  and  f  if  going  from  two  to  six  terms  for  the  same 
model  parameters  are  shown  in  Fig.  4.  The  dashed  lines  are  the  two  term  results. 
The  two  term  results  underestimate  the  final  values  at  both  the  low  and  high  energy 
ends  while  overestimating  at  the  intermediate  energies.  This  leads  to  the 
conclusion  that  the  reaction  rate  coefficients  calculated  using  the  two  term 
expansion  tend  to  underestimate  their  true  values. 
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W  (cm/sec) 


Figure  2 


The  transport  coefficients  W  and  DT/p  as  a  function  of  E/N  calculated 

—16  2 

for  Np=2  and  Np=6  in  the  model  atom  with  k=10  x  10  cm  and  e.=,2eV. 


18 


Dt//a  (eV) 


ENERGY  (eV) 


10 


ENERGY  (eV) 

Figure  4.  Comparison  of  f  and  f.  from  the  two  term  and  six  term  solutions. 

O  1 

The  solid  lines  are  the  Np“6  solutions  and  the  dashed  lines  are  for  Np*2. 
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Methane 


The  cross  sections  adopted  for  these  CH^  calculations  are  as  given  by 

g 

Kleban  and  Davis*  These  cross  sections  shown  as  Fig.  5  are  only 

meant  to  represent  the  gross  features  of  the  transport  data  and  were  not  rigorously 
derived.  We  choose  CH^  for  a  detailed  analysis  because  of  the  large  value  of 
the  Inelastic  cross  section  in  the  region  of  the  Rarosauer  minimum  in  the  elastic 
cross  section.  Based  on  the  previously  assumed  criterion  of  a  small  ratio  of  inelastic 
to  elastic  cross  sections,  the  two  term  analysis  using  these  cross  sections  may  be 
expected  to  exhibit  severe  departures  from  the  multi-term  analysis.  This  has 
indeed  been  the  case. 

The  individual  expansion  coefficients  f^  for  i*=0  to  5  are  shown  in  Fig.  6 

for  E/N  *=  2.42  Td.  The  main  features  seen  in  the  case  of  the  model  atom  can 

also  be  seen  here,  i.e.,  the  zeros  and  two  maxima  in  the  higher  order  coefficients. 

Ficrure  7  shows  f  _  and  f  -  for  N  -2  and  N  —f\ .  As  for  rhp  modol  prom  th&  rwo  farm 
u  x  r  r 

solution  tends  to  overestimate  fQ  and  f^  for  the  intermediate  energies  and 
underestimates  their  values  at  higher  energies. 

Figure  8  shows  results  of  calculations  of  W  and  DT/y  for  N=2  and  N=6 
over  a  range  of  E/N  around  the  maximum  in  W.  In  no  case  is  the  two  term  value 
of  D^/y  within  5%  of  the  multi-term  results.  Again  as  for  the  model  atom, 
the  change  in  going  from  Np-2  to  4  is  larger  than  that  from  Np=4  to  6  or  8. 

We  conclude  that  the  two  term  analysis  is  inaccurate  for  CH^  at  this  range 
of  E/N  and,  although  the  error  seems  to  be  decreasing  slightly  as  E/N 
is  increased,  there  is  probably  no  region  where  the  two  term  solution  in  CH^ 
is  valid  to  better  than  5%* 

Monte  Carlo  calculations  we  have  made  in  CH^  verify  our  multi-term  results 

9 

for  the  transport  coefficients*  Lin,  et  al*  have  developed  a  moment  method 
for  solving  the  Boltzmann  equation  which  they  applied  to  CH^  using  the  present 
cross  sections*  We  compare  well  with  their  published  values  at  3*9  Td* 
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Figure  6.  The  first  six  Legendre  expansion  coefficients,  f^  for  i*0  to  5, 

as  a  function  of  energy  for  methane  at  2.42  Td. 
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Figure  7..  Comparison  of  fQ  and  f ^  from  the  two  term  and  alx  term  solutions. 
The  solid  lines  are  Np  ■(  solutions  and  the  dashed  lines  are  for  Np*2. 


Figure  8.  The  transport  coefficients  W  and  Dt/(j  as  a  function  of  E/M  calculated 


for  Npa2  and  Np*6  in  methane. 
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SECTION  VI 


NITROGEN  RESULTS 


The  multi-term  Boltzmann  code  has  been  used  to  calculate  transport  parameters 

in  N2  at  1,  40,  and  100  Td.  As  input  to  the  program,  we  use  the  set  of  cross 

10 

sections  in  derived  by  Phelps,  Our  point  here  is  to  evaluate  the  error 
introduced  in  the  calculated  W  and  D,pN  due  to  the  two  term  expansion  approximation 
and  to  compare  with  other  methods. 

We  have  made  the  calculations  in  using  three  different  approximations 
in  the  Boltzmann  code  and  using  a  Monte  Carlo  method.  These  different  approxi¬ 
mations  were  used  with  the  Boltzmann  approach  to  get  some  feeling  for  where  each 


was  more  appropriate  and  easier  to  implement.  These  three  approximations 
were  all  possible  using  our  multi- term  code.  We  have  already  discussed  two 
of  these  in  some  detail,  the  two  term  and  the  multi-term  approximations. 

A  different  approximation  Is  due  to  Baraff.^ 

H 

The  Maximum  Anisotropy  Truncation  (MAT)  scheme  proposed  by  Baraff  begins  by 
assuming  that  the  distribution  function  is  highly  elongated  and  an  expansion  in 
Legendre  functions  is  very  slowly  converging  if  at  all.  If  the  angular  dependence 
of  f(v)  can  be  represented  by  a  delta  function,  the  relation  between  adjacent 


Legendre  expansion  coefficients  is 


f4<E>  •  fi-i<E) 


Rather  than  assuming  that  the  higher  order  expansion  coefficients  can  be 
neglected  (f^-0  for  i>l)  as  in  the  traditional  two-term  Boltzmann  analysis,  the 
MAT  truncation  or  uncoupling  is  based  on  Eq.  (21).  Baraff  solves  the  two-term 
Boltzmann  equation  using  Eq.  21  to  substitute  for  f2  ln  terms  of  f^  in  Eq.  (6). 
We  have  implemented  the  MAT  in  our  code  and  have  solved  the  Boltzmann  equation 


keeping  two  terms  and  then  six.  If  the  higher  order  coefficients  are  Indeed 

negligible.  It  should  make  no  difference  how  the  series  Is  truncated  but  one 

method  may  require  less  terms  than  the  other  to  converge.  As  will  be  seen, 

the  series  truncation  by  setting  f2*0  does  yield  two-term  values  closer  to 

the  converged  results  than  the  MAT  scheme. 

The  Monte  Carlo  calculations  were  performed  using  the  code  provided  by 

Ivan  Reid  and  developed  at  the  Australian  National  University.^ 

Table  I  shows  the  results  of  our  calculations  of  W  and  D^N  using  the  four 

approaches  discussed  above.  The  first  three  rows  of  values  for  W  and  DTN  were 

computed  by  setting  f^  “  0  in  equation  6.  A  check  of  the  Np=2  results 

*  3 

were  made  using  the  backward  prolongation  method.  As  in  the  previously  discussed 
model  cases,  the  change  in  the  values  of  W  and  D^N  in  going  from  Np*2  to 
Np*=4  are  much  greater  than  from  Np~4  to  Np*6.  The  differences  between 
the  Np*2  and  Np=6  values  are  greater  for  higher  values  of  E/N,  as  expected, 
and  are  greater  for  D^N  than  for  W. 

The  fourth  and  fifth  rows  are  values  calculated  using  the  MAT  assumption. 

There  is  a  large  difference  between  the  two  truncation  methods  for  Np=2 
and  a  very  small  difference  for  the  Np^6  cases.  This  is  to  be  expected  if 
the  series  is  converging.  The  MAT  Np*2  values  are  in  substantial  disagreement 
with  the  MAT  Np*6  values. 

The  last  row  contains  values  of  W  and  D^N  calculated  from  a  Monte  Carlo  method 

There  are  some  discrepancies  between  the  Boltzmann  results  and  the  Monte  Carlo 
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results  at  an  E/N  of  4.0*10  Vcm  .  We  have  reason  to  suspect  that  the  Monte  Car] 
code  may  be  in  error  at  the  higher  E/N  values.  For  one  thing,  between  collisions 
in  the  Monte  Carlo  code,  the  electrons  travel  in  straight  lines.  The  effect 


TABLE  1 


N^  Drift  Velocity  -  W 


E/N 

1  x  10"17  V  cm2 

4  x  l(f 16  V  cm2 

1 x  10"15  V  cm2 

Hp-2 

4.118x10^cn  sec  1 

5.640  x 106cm  sec"1 

1.111  *  107cra  sec'1 

V4 

4.118 

5.560 

1.096 

Np-6 

4.118 

5.558 

1.096 

MAT  Np-2 

4.149 

6.177 

1.286 

MAT  Np“6 

4.118 

5.559 

1.105 

Monte  Carlo 

4.12 

5.42 

1.083 

Nj  Diffusion  Coefficients  -  D^,N 


E/N 

1  X  10'17  v  cm2 

4  x 10" 16  V  cm2 

1 x  10"15  V  cm2 

V2 

1. 323 x lO^^cm  ^sec  ^ 

1*744  *  lO^cm  ^sec  ^ 

2*035  x  lO^cm  ^sec  ^ 

Np«4 

1.317 

1.624 

1.822 

Np«6 

1.319 

1.630 

1.850 

MAT  Np*2 

1.326 

1.758 

1.953 

MAT  Np*6 

1.317 

1.612 

1.852 

Monte  Carlo 

1.305 

1.604 

1.901 
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of  the  field  is  to  curve  the  electron  trajectories.  The  straight  line  paths 
nay  tend  to  reduce  the  calculated  drift  velocities ,  and  this  is  what  is  seen 
when  comparing  with  the  Boltzmann  code. 


At  this  point  it  is  not  too  instructive  to  compare  computing  speeds  among 
the  various  methods.  First,  very  little  effort  has  been  made  to  optimize  the 
speed  requirements  in  the  Boltzmann  code.  The  Monte  Carlo  calculations  are 
very  time  consuming,  but  an  exact  figure  on  the  savings  using  the  Boltzmann  code 

is  not  available  yet.  Presently  factors  of  10  to  1000  in  speed  are  saved  by 
using  the  Boltzmann  code.  Another  point  making  a  comparison  difficult  is  the 
fact  that  we  are  using  the  vector  processing  facilities  of  the  CRAY-1  computer. 


It  is  hard  to  estimate  what  execution  times  would  be  without  that  option. 


SECTION  VII 


CONCLUSIONS 

The  determination  of  low  energy  electron  scattering  cross  sections  from  swarm 
experiments  has  relied  almost  exclusively  on  the  two  term  Boltzmann  analysis 
up  to  now  for  the  calculation  of  W  and  DT/u.  This  has  generally  been  assumed 
valid  for  cases  where  the  inelastic  cross  sections  are  much  smaller  than  the 
elastic.  We  have  presented  here  a  method  for  generalizing  the  Boltzmann  analysis 
to  include  an  arbitrary  number  of  terms  in  the  Legendre  expansion  of  the  distri¬ 
bution  function.  The  method  has  been  applied  to  two  model  cases  and  to  nitrogen 
as  an  example  of  a  real  gas  with  a  large  number  of  inelastic  processes. 

We  find  that,  in  general,  the  error  introduced  in  the  transport  parameters 
calculated  using  the  two  term  approximation  is  greater  at  higher  ratios  of  inelastic 
to  elastic  cross  sections.  In  the  case  of  methane,  the  error  decreases  with 

decreasing  ratio  of  cross  sections  yet  increasing  E/N.  The  error  in  DTN  for 

any  particular  case  is  always  much  larger  than  for  W. 

Agreement  between  other  methods  of  our  multi-term  results  are  good. 

Comparison  between  the  backward  prolongation  two  term  calculation  and  our  code 

with  Np=2  is  exact  for  all  cases  investigated.  Our  converged  values  of  transport 

coefficients  compare  well  with  Monte  Carlo  results  within  quoted  accuracies 

except  for  the  higher  E/N  cases  in  The  reason  for  the  discrepancies  is  not 

understood.  Finally,  since  this  project  began  another  multi-term  Boltzmann  method 
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based  on  moment  techniques  has  appeared  in  the  literature.  We  compare  our  results 

in  CH^  at  3.9  Td  with  the  values  published  for  W  and  DTN  as  a  function  of  Np 
and  find  excellent  agreement. 
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We  have  been  concerned  here  that  part  of  the  swarm  experiment  analysis 

in  which  the  Boltzmann  equation  is  solved  for  the  electron  energy  distribution 

function  and  hence  the  transport  coefficients.  The  iterative  determination 

of  cross  sections  by  comparison  of  calculated  and  measured  values  of  W  and  D/p 

Is  yet  another  issue  about  which  only  general  remarks  may  be  made.  Kumar, 

12 

Skullerud,  and  Robson  note  that  in  order  to  glean  maximum  information 

from  the  potentially  extremely  accurate  swarm  measurements,  we  should  aim 

for  accuracies  of  .1%  and  1%  for  W  and  respectively,  in  the  calculations. 

We  have  shown  here  that  such  accuracies  can  not  be  achieved  using  the 
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two-term  approximation.  It  is  possible  to  show  that  in  helium,  a  2% 
variation  in  the  momentum  transfer  cross  section  leads  to  about  a  2%  change 
in  the  calculated  drift  velocity  over  a  certain  range  of  E/N,  but  it 
should  be  noted  that  the  correspondence  between  change  in  cross  sections 
and  change  in  drift  velocity  is  certainly  not  in  general  true. 

The  validity  of  the  two  term  approximation  has  been  debated  and  questioned 
for  some  time.  If  accuracies  of  .1%  and  1%  or  less  in  W  and  D^N,  respectively, 
are  desired  from  the  calculations,  the  two  term  approximation  is  of  limited 
validity.  In  N2  at  40  Td,  the  two  term  W  and  D^N  values  differ  from 
the  six  term  values  by  1.4  percent  and  6.9  percent,  respectively.  At  100  Td, 
the  differences  are  1.6  percent  and  10  percent,  respectively.  However,  for  1  Td, 
the  two  term  values  seem  quite  accurate.  In  CH^  the  two  term  DTN  is  greater 
than  5  percent  different  from  the  six  term  value  for  all  values  of  E/N 
Investigated. 

Attempts  to  quantify  the  error  introduced  by  the  two  term  approximation 
have  met  with  limited  success.  We  plan  to  continue  these  efforts  in  this  direction, 
as  some  criteria  for  when  the  two  term  approximation  is  valid  would  be 
very  useful. 
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A  variety  of  problems  may  now  be  addressed  using  the  code.  One  interesting 
problem  is  the  behavior  of  the  high  energy  tail  of  the  distribution  function 
as  a  function  of  the  number  of  Legendre  functions  in  the  solution.  This  will 
hopefully  resolve  some  of  the  problems  encountered  in  determining  electronic 
excitation  or  ionization  rate  coefficients  simultaneously  with  the  transport 
coefficients.  Finally,  with  the  method  presented  here,  we  may  now  return  to 
the  original  problem  of  the  analysis  of  the  swarm  experiments  in  water  vapor 
to  obtain  more  accurate  values  of  the  low  energy  electron-water  molecule 
scattering  cross  sections. 


32 


REFERENCES 


1.  J.  R.  Rumble ,  Jr.,  E.  C.  Beaty,  and  L.  C.  Pitchford,  In  Proceedings  of  the 

Sixth  International  Codata  Conference  (Sicily,  May  1978),  Bertrand 
Dreyfus,  editor,  Permanom  Press,  Oxford,  243  (1979). 

2.  H.  R.  Skullerud,  J.  Phys.  B  2,  696  (1969). 

3.  L.  S.  Frost  and  A.  V.  Phelps,  Phys.  Rev.  127,  1621  (1962). 

4.  N.  L.  Schryer ,  Bell  Laboratories  Computing  Report  Number  52  (1976). 

5.  C.  de  Boor,  J.  Approx.  Theory  <>,  50  (1972). 

6.  C.  de  Boor,  in  "Spline  Functions  and  Approximation  Theory,"  A.  Meir  and 

A.  Sharma,  ed.,  ISNM  vol.  21,  Birkhauser  Verlag,  Basel  (1973),  57. 

7.  Ivan  Reid,  thesis,  Australian  National  University,  1978. 

8.  P.  Kleban  and  H.  T.  Davis,  Phys.  Rev.  Lett.  3^,  456  (1977). 

O  C  T  T  ft  P  T>_l* —  1  n  »  ««  *  — -  -•  *%  »  ^  V 

•  *  -  — •  - »  wttu  •  n.  iMOUU|  O  •  V/UClIl  •  r  liyb*  /  A  »  J  )  • 

10.  A.  V.  Phelps,  D.  Levron,  and  K.  Tachibana,  in  Proceedings  of  the  Xlth  Inter¬ 

national  Conference  of  the  Physics  of  Electronic  and  Atomic  Collisions 

(Kyoto,  1979),  K.  Takayanagi  and  N.  Oda,  editors. 

11.  G.  A.  Baraff,  Phys.  Rev.  133,  A26  (1964). 

12.  K.  Kumar,  H.  R.  Skullerud,  and  R.  E.  Robson,  "Kinetic  Theory  of  Charged  Particle 

Swarms  in  Neutral  Gases,"  to  be  published. 

13.  B.  Bederson  and  L.  J.  Keiffer,  Rev.  Mod.  Phys.  43,  601  (1971). 

THIS  PAGE  IS  BEST  QUALITY  PRACTICASll 

F 't  1  v.'*-  *  ..  a;  juC 


33 

eU.S.Oov«rnm*nt  Printing  Offlcti  1980  —  657-084/723 


